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We have investigated the structural sequence of the high-pressure phases of silicon and germanium. 
We have focussed on the cd — > /3-tin — > Imma — > sh phase transitions. We have used the plane-wave 
pseudopotential approach to the density-functional theory implemented within the Vienna ab-initio 
simulation package (VASP). We have determined the equilibrium properties of each structure and 
the values of the critical parameters including a hysteresis effect at the phase transitions. The 
order of the phase transitions has been obtained alternatively from the pressure dependence of the 
enthalpy and of the internal structure parameters. The commonly used tangent construction is 
shown to be very unreliable. Our calculations identify a first-order phase transition from the cd to 
the /3-tin and from the Imma to the sh phase, and they indicate the possibility of a second-order 
phase-transition from the /3-tin to the Imma phase. Finally, we have derived the enthalpy barriers 
between the phases. 

PACS numbers: 61.50.Ks, 64.70.Kb, 71.15.Nc, 71.20.Mq, 81.30.Dz, 81.40.Vw 



I. INTRODUCTION 

Although the structural sequence of the high-pressure 
phase transitions in Si and Ge has been considered to be 
well known since two decades, a new high-pressure phase 
was found experimentally between the body-centered 
tetragonal structure with two basis atoms (which will be 
referred to the following as /3-tin or BCT phase) and the 
simple-hexagonal (sh) structure in Sii*£ and Gei^ After 
its space group this intermediate phase was called Imma 
phase and corresponds to a body-centered orthorhombic 
Bravais lattice with two atoms in the unit cell (which 
will be called in the following BCO). Thus, the exper- 
imental structural sequence is found to be cd — > /3-tin 
— > Imma — > sh — > . . ., where cd indicates the cubic- 
diamond structure. A BCO structure was first theo- 
retically proposed for Si by Needs and Martini about 
twenty years ago, whereas Lewis and Cohen£ predicted 
exactly the Imma phase for Ge a few years before the ex- 
perimental confirmation. The first theoretical investiga- 
tions for Si and Ge indicated for both, the /3-tin— >Imma 
and the 7mma->sh continuous second-order phase 

transition*^ At variance, the experiments show first- 
order phase transitions with a discontinuity of the vol- 
ume in both substances and a hysteresis effect in Sii£ For 
Ge the first-order transition was detected only for the 
/3-tin— >Imma transition, 3 because the Imma— >sh transi- 
tion was expected to be beyond the experimentally ac- 
cessible pressure range. The transition pressures in Ge 
are expected to be higher than the corresponding ones 
in Si, because of the strong repulsive character of the 
atomic potential in Ge due to the presence of d elec- 
trons in the core. 8 From the theoretical point of view, re- 
cent total-energy calculations show two first-order phase 
transitions /3-tin— *Imma and Imma— >sh both for SiS and 
GeiiS A second-order phase transition was found for the 
ft-tin—* Imma case from elastic instabilities^ 

First information on the order of the phase transition 
can be given by group-theoretical arguments*^ If the or- 



der of the point group of one phase is one half of the or- 
der of the point group of the other phase, then the phase 
transition can be of second order. If the point group of 
one phase is one third of the order of the point group 
of the other phase, then the phase transition must be of 
first order. The order of the point group m3m of the cd 
structure (space group 7 h , Fd3m) is 48. Since the or- 
der of the point group A/ mmm of the /3-tin phase (space 
group D^t, IAi/amd ) is 16, the cd^ /3-tin transition 
has to be of first order. The order of the point group 
mmm of the Imma phase (space group D^, Imma) is 8, 
thus the /3-tin— >7mma transition can be a second-order 
one. Furthermore, the order of the point group 6 /mmm 
of the sh structure (space group D\ h , PQ/mmm) is 24, 
the sh— > Imma and also the Imma— >sh transitions are of 
first order. Following the group-theoretical analysis none 
of the available total-energy calculations have described 
correctly the order for the phase transitions of the whole 
sequence, neither for Si nor for Ge. 

This paper is organised as follows: In Section [D] we 
give a short overview of the methods on which our calcu- 
lations are based. In Section ITTll we report the results for 
the equilibrium parameters of these phases and compare 
them with the experimental volume dependence. Tran- 
sition pressure ranges are presented in Section llVI The 
order of the phase transitions is analysed in Section 
and the enthalpy barriers between the phases are deter- 
mined in Section IVll All results are summed up and dis- 
cussed in Section IVIII The conclusions are summarised 
in Section 



II. METHOD 

Our calculations have been performed using the plane- 
wave pseudopotential approach to the density-functional 
theory (DFT)i2*i^ implemented within the Vienna ab- 
initio simulation package (VASP)^ For the exchange- 
correlation potential we have used the generalised- 
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gradient approximation (GGA) by Perdew and Wang 16 
for Si and Ge, and, because of poorer results, the local- 
density approximation (LDA)1Z^ for Si only. The in- 
teraction with the ion cores has been described by ul- 
trasoft pseudopotcntialsii2i2S The wave functions have 
been expanded in terms of plane waves up to a kinetic- 
energy cutoff of 270 eV (410 eV) for Si (Ge). This choice 
provided an error smaller than 0.5 kbar (0.2 kbar) for 
Si (Ge) to the pressure according to the Pulay stress. 21 
The Pulay stress arises from using an incomplete plane- 
wave basis set and causes errors in the stress tensor 
and, correspondingly, in the pressure. Convergence also 
required a 18x18x18 (24x24x24) mesh of Monkhorst- 
Pack points 22 which yielded 864 (1962) k-points in the 
irreducible wedge of the Brillouin zone for Si (Ge). All 
high- pressure phases are metallic. Thus, we have used 
a Methfessel-Paxton smearing 2 ^ with a width of 0.2 eV. 
In order to determine the total energy for a given vol- 
ume or structure we have minimised the energy with 
a conjugate-gradient algorithm. 24 - 25 - 26 The equilibrium 
properties have been obtained by fitting the values of the 
energy as a function of the volume alternatively to the 
Vine1*2£ and to the Murnaghan 2 ^ equation of state. 



III. STRUCTURAL RELAXATION AND 
EQUILIBRIUM PROPERTIES 

In this section we present our results for the volume de- 
pendence of the structural parameters and compare them 
with the experimental ones. The structural parameters 
are more useful indicators of phase transitions than the 
often-used tangent construction to the equation-of-state 
curves, the difficulties of which will be illuminated. 

We have focussed our attention to the structural se- 
quence /3-tin— >Imma— >sh. The calculations for the cd 
phase have been performed just for the sake of compar- 
ison. The most general structure which is compatible 
with all four phases is the body-centered orthorhombic 
cell with a basis of two atoms at (0,0,0) and (0,0. 56, Ac) 
(BCO structure). The other structures are special cases: 
For the BCT structure (/3-tin) one has the lattice con- 
stants a = b and the internal parameter A = 0.25; for 
the sh structure one has b = V3c and A = 0.5; for the 
cd phase one has b = a, c — V2a, and A = 0.25. The 
pure Imma phase is realised as a BCO structure with 
a b =/= c =/= a and 0.25 < A < 0.5. The minimization of 
the total energy has been performed using the following 
procedure: For all phases we have started from the BCO 
structure and minimised the energy as a function of the 
structural parameters with the constraints appropriate 
to the various phases. In this way we have avoided an 
energy offset which arises from using different cells. For 
example, in the case of the cd structure we have imposed, 
according to the previous considerations, b = a, c = y/2a, 
and A — 0.25, and then we have minimised the total en- 
ergy as a function of a. In the special case of the sh struc- 
ture, we could keep only the value of A fixed at 0.5, and 



we have relaxed the value of b/c (and also c/a) by min- 
imization of the energy; the structure which is obtained 
after relaxation will be indicated in the following as SH. 
Indeed, converged calculations should be yield b/c = \/3 
(SH = sh) after relaxation. For the Imma phase, all lat- 
tice and internal parameters of the BCO structure were 
allowed to relax. 
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FIG. 1: Volume dependence of the structural parameters of 
the BCT, BCO, and SH structures for Si within GGA (upper 
panels) and LDA (lower panels). The dashed line indicates 
the ideal v alue of b/c = V3. The experimental data are taken 
from Refs. 1 2 3CM3lll32 



The relaxed parameters as a function of the volume are 
displayed in Figs-^andElfor Si and Ge, respectively. The 
results from calculations performed withing the GGA and 
the LDA are quite similar for Si. Since it turned out, 
that the pressures calculated by GGA are more accurate 
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in Si, only GGA calculations have been performed for 
Ge. As can be seen in Fig. the relaxation of b/c for 
the SH phase of Si does not match the ideal value. At 
variance, the calculation for Ge does reproduce this value. 
This indicates that for a full convergence a higher kinetic- 
energy cutoff and/or a larger grid of special points should 
have been used for Si. However, the condition b/c = \/3 
is almost fulfilled in the region of stability of the sh phase. 
Furthermore, the energy obtained from a calculation in 
the same cell with b = V3c fixed and c/a relaxed differs 
from the value obtained using the previous procedure by 
less than 0.7 meV. For a more detailed discussion see 
Ref.@ 

The agreement between our data and those experimen- 
tal ones which are not included in the figure is good: 
The reported value of c/a = 0.552 for the /3-tin phase 
of SiS is well within the range of our values. The ex- 
perimental values of b/a = 0.917 and c/a = 0.530 at 
the volume V = 24.43 A 3 of Ref. Ei agree with our 
data. Furthermore, the measured value of 24.22 A 3 for 
the volume of the sh phase of Si 3 ^ is also in the stabil- 
ity range we have found for this phase. For Ge in the 
BCT structure, a c/a ratio of 0.548-0.554 was detected 
experimentally 32 ' 36 ! 3 ^ 3 ^ 3 ^ corresponding to a volume re- 
gion around 33 A 3 ; in this region, our c/a ratio is 0.546- 
0.549. For the sh phase of Ge, a b/a ratio of 0.93 was 
measured 3 ^ which also agrees with our data. 



0.50 



,.: L. \ ■ . ■. .■ . — , .■ — . . . '-. ,'. ,' .. — w~~ — ,- . ,'\ / . ,"— 



<] 



0.95 

0.90 
1 .90 

1 .80 

1 .70 



Ge (GGA) % t 



• EXP □ BCT 

o BCO a SH 



0.25 LrinnnnmnnnnmnnnnmnnnnmnnnnmnnnRQariftnr^- 

0.55 

^□□□□□□□□□□□□□□□□□ nnnnnnnnnn 8 B B DBan 

' 4 !- *oo°o 

\AA^^6QS 8 2 AAAAAAAAAAAAA 

0.53 

j qq ^uuuuuuuuuuuuuuuuuuumuuuuwuuujjtyyuuuiy- 



,QO> 



)CQQ0(SQfflfi BA AAAAAAAAAAAAAAAAAAAAAAA 




23.0 24.0 25.0 

Volume (A 3 ) 



26.0 



FIG. 2: The same of Fig. □ for Ge within GGA. The experi- 
mental data are taken from Ref. 0- 

Figures ^ an d 12] show that the relaxation of the gen- 
eral BCO structure leads to the BCT structure for large 
volumes and to the SH structure for small ones. In the in- 
termediate region, extending from about 26.5 to 28.0 A 3 
for Si and 23.5 from 25.5 A 3 for Ge, we localise the sta- 
bility range of the pure Imma phase. Whereas the transi- 



tion between the Imma and the sh phase is characterised 
by a discontinuity in most of the structural parameters, 
the Imma and the /3-tin phase are quite difficult to dis- 
tinguish by an inspection of Figs. ^ and [21 A clearer 
picture is obtained from the electronic structure and the 
Fermi energy (see Appendix lAl for details). From this 
procedure, the stability range of the Imma phase is from 
26.6 to 28.3 A 3 for Si (GGA), from 26.6 to 28.0 A 3 for 
Si (LDA), and from 23.7 to 25.4 A 3 for Ge. 

The standard procedure for calculating the volume 
change and the pressure at the transition is the so-called 
common-tangent construction. It relies on the fact that, 
in the case of a first-order phase transition, two struc- 
tures with the same pressure p and enthalpy H but dif- 
ferent total energy E(V) coexist. The volume derivative 
of the total energy is the negative pressure. Therefore, 
the slope of the common tangent to the two energy curves 
gives the pressure for the transition between the corre- 
sponding structures. Due to the fact that calculations are 
performed at a finite number of volumes, a fit of an ana- 
lytical expression to a continuous energy curve is needed 
in order to obtain a continuous derivative. 



-9.80 



-9.90 - 



=S -10.00 



-10.10 



-10.20 



, Si (GGA) 




• Imma 






~ BCO 




\ 


° BCT 






a SH 
















Spa AAAAAAAAAAAAA. 







26 



27 



28 29 
Volume (A 3 ) 



30 



31 



FIG. 3: Energy vs. volume for the involved structures of Si 
within GGA (see text for explanation). 

In Fig. |3|we show the energy of the relaxed structures 
for Si as a function of volume. The volume region where 
the pure Imma phase is stable does not include any min- 
imum of the total energy E. This makes a fit to any 
equation of state for the Imma phase slightly unstable, 
which results in an unreliable determination of the tran- 
sition pressures and of the volume change through the 
common-tangent construction in this region. The results 
from a fit of the VineliSi and Murnaghan 28 equation of 
state are shown in Table [i] for Ge as well as a comparison 
of the LDA and GGA calculations for Si. 

Typically, the GGA tends to overestimate the equilib- 
rium volume whereas the LDA tends to underestimate it. 
For the cd phase of Si the experimental equilibrium vol- 
ume at room temperature i^^^L^^ 40.05 A 3 , which 
lies between the LDA and GGA ones, with a deviation of 
about 2%. The experimental value for cd Ge [^J^JUk 
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TABLE I: Equilibrium parameter determined by the Mur- 
naghan (Mur) and Vinet (Vin) equation of state for the cd, 
/3-tin, Imma, and sh phase of Si and Ge. The equilibrium 
volume Vb is in units of A 3 and the bulk modulus Bo in kbar. 
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/3-tin 
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Bo 
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B'o 


4.05 


4.61 


4.20 


4.52 


3.50 


5.24 



between 45.00 and 45.31 A 3 , with a deviation from our 
result of less than 6%. Because the other phases exist 
only under pressure, experimental values do not exist to 
which our equilibrium parameters could be compared. 

If the standard common-tangent construction is used 
in order to extract the transition parameters, an accurate 
determination of both, equilibrium volumes and curva- 
ture of the energy vs. volume curve, is needed. This 
applies even more so in the present case in which E(V) 
curves and equilibrium volumes of the /3-tin, Imma and sh 
phases are very close together. Even though the present 
calculation are numerically involved the resulting data of 
Tableware not precise enough. For example, in all calcu- 
lations (LDA and GGA) for Si the fitted (Vinet and Mur- 
naghan) energy curves of the Imma phase fall below the 
ones of the /3-tin phase, even for large volumes. In these 
cases the Imma phase is more stable than the /3-tin phase. 
Thus, within this construction cither a phase transition 
from the cd to the /3-tin phase can be observed, if one 
neglects the Imma phase, or a phase transition from the 
cd to the Imma phase can be observed without a locally 
stable /3-tin phase. A direct transition from the /3-tin to 
the Imma phase cannot be detected. Similar drawbacks 
are found for Ge. Furthermore, the method is applicable 
only to the investigation of first-order phase transitions. 
For all these reasons, it is very difficult to extract the 
reliable information on the transitions using the tangent 
construction. In the following section, we consider an al- 
ternative procedure for analysing the transitions which 
does not contain the drawbacks of the above one. 



IV. DETERMINATION OF THE TRANSITION 
PRESSURES AND HYSTERESIS EFFECT 

In this section we present results for the pressure de- 
pendence of the structural parameters in order to give 
an estimate of the transition pressures by considering the 
stability range of the Imma phase. 

The pressure p in a crystal can be directly calculated 
from the stress tensor <x. Within the BCO structure 
the off-diagonal components of the stress tensor vanish. 
Therefore, the most general form of the stress tensor of 
our structures is 

(p. \ 

° = - Py \ ■ (1) 

V p. J 

The results of the previous section correspond to hydro- 
static pressure 

P = Px=Py=Pz , (2) 

and, in our results of Section llTTl Eq. J5J is fulfilled within 
an error of 0.1 kbar for Si and 0.05 kbar for Ge, which 
is less than the pressure error due to the Pulay stress. 
For each of the relaxed BCO structures as described in 
Section ITTT1 we have calculated the hydrostatic pressure 
corresponding to the structural parameters as shown in 
Figs. ^ an d HI The pressure dependence of the struc- 
tural parameters is presented in Figs. 01 and |SJ There 
is a clear discontinuity of the parameters of the Imma 
phase and the sh phase. The values of the parameters of 
the pure Imma phase are highlighted as special values of 
the relaxed BCO structure. For both Si and Ge, the /3- 
t'm^Imma phase transition seems to be continuous. For 
Si, the relaxed parameters calculated within LDA and 
GGA are very similar, but the corresponding pressures 
differ as shown in Fig. 01 In contrast, the volume range of 
stability in both cases is similar (see Fig.^l. Usually, the 
GGA transition pressures are closer to the experimental 
ones than the LDA values. Therefore, we will focus in 
the following mainly on the GGA results. 

The assignment of the phases from the data points in 
Figs. 0] and [3] are taken from Section IIIII In the case of 
the GGA calculation for Si we have found the following 
sequence with pressure increase: For the /3-tin— > Imma 
phase transition, the transition pressure is approximately 
between 104 (last relaxation of the system to a BCT 
structure) and 109 kbar (first relaxation of the system 
to a BCO structure with A ^ 0.25). The next transi- 
tion is the one from the Imma to the sh phase. As can 
be seen in Figs. 0] and [SI there is an overlap of different 
structures for certain values of the pressure, which can 
be interpreted as a hysteresis effect. That means that the 
transition value for increasing and decreasing pressure is 
different. This effect was first observed first in the ex- 
periment for Sii£ We estimate the hysteresis effect for Si 
as follows: The last appearance of the pure Imma phase 
(A ^ 0.5) without any SH structure at the same pressure 
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TABLE II: Summary of the estimated transition pressures 
at increasing (f) and decreasing (J.) pressure for the phase 
transitions between the /3-tin, Imma and sh phase. 
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FIG. 4: Dependence of the structural parameters A and b/a 
on the hydrostatic pressure for the relaxed phases for Si with 
GGA (full and empty diamonds) and LDA (full and empty 
circles). The range of the pressure hysteresis of the overlap- 
ping phase regions are marked by vertical lines. 
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FIG. 5: The same of Fig. H but for Ge for GGA only. 



is at 186 kbar, the first appearance of the sh phase is at 
191 kbar (downstroke) , the last appearance of the pure 
Imma phase is at 197 kbar and the first appearance of 
the sh structure without any parallel Imma phase is at 
198 kbar (upstroke). The transition pressures estimated 
in this way are listed in Table ILT1 

Since we have calculated the lowest energy structure, 
the values in Table [H] give an estimate of the hysteresis 
effect. An alternative approach to the hysteresis effect 
will be described in Section IVll 



The thermodynamical potential for a pressure-induced 
phase transition is the Gibbs free energy. Neglecting tem- 
perature effects, the Gibbs free energy is equal to the 
enthalpy H = E + pV. For calculating the enthalpy, we 
have used the energy E(V) of the relaxed structures with 
a volume V and the hydrostatic pressure p from Eq. 

At a given pressure, the stable structure is the one with 
the lowest enthalpy. Therefore, the transition pressure pt 
for a transition A-^B is defined by H A (p t ) = H B (p t ). 
For a first-order phase transition the enthalpy curves 
H A (p) and H B {p) of the two phases cross at p t , and their 
derivatives AH /dp at p t are different for the two phases. 
For a second-order transition, the second derivatives of 
the enthalpy curves at pt for the two phases are different, 
whereas the first derivatives are equal. In this case the 
corresponding enthalpy curves have a boundary point at 
Pt- 

In the following, we analyse the transitions in order to 
extract the corresponding transition pressure and order. 
The enthalpy H(p) and its first derivative are presented 
in Fig.|S|for Ge. The enthalpy curves have been reduced 
by the value of the BCT structure for the /3-tin— ► Imma 
transition, 



AH(p) = H(p)-H BCT (p) 

and by the value of the SH structure for the Imma- 
transition, 

AH(p)=H(p)-H sii (p) . 



(3) 
->sh 



(4) 



The derivative indicated by black dots in the lower panels 
of Fig. has been calculated by finite differences, the 
straight solid line is a linear fit of these data, obtained 
by using only points belonging to the numerically stable 
regime (790-1030 kbar). An integration of this linear 
function is fully consistent with the enthalpy values in the 
upper panels. The enthalpy differences of the points in 
the pressure range p < 790 kbar stem from the numerical 
noise of our calculation. 

For the /mma-*sli transition (right panels of Fig.[5J) the 
two enthalpy curves cross at the transition pressure pt = 
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FIG. 6: Enthalpy vs. pressure (upper panels) and derivative 
of the enthalpy vs. pressure (lower panels) for Ge reduced 
by the corresponding curves of the BCT (right panels) and 
the SH (left panels) structure. The solid lines are quadratic 
(upper panels) and linear (lower panels) fits. 



1009 kbar. Therefore, this phase transition is clearly of 
first order. The change of the derivatives at the transition 
point is the volume change AV(pt) and can be directly 
read from the lower panel of Fig. El 

For /3-tin— > Imma (left panels of Fig. the enthalpy 
curves are very close to each other in the region where 
the phase transition is expected. From this point of view, 
transitions of both first and second order are conceivable. 
We have thus analysed our data assuming alternatively 
a first- and a second-order transition. 

Assuming a second-order transition the transition pres- 
sure is obtained by looking at the crossing of the dH/dp 
curves. 

A precise value of the /3-tin— >Imma transition pressure 
assuming first order can not be determined by a crossing 
of the enthalpy curves due to the numerical accuracy of 
our calculation. We are only able to give a rough estimate 
of the pressure at which the transition can occur due to 
an extrapolation of the more reliable data points away 
from the transition pressure. The results for Si have been 
obtained in an analogous way. 

For both Si and Ge, the values of pt and of the volume 
changes between the low- and the high-pressure phase 
are shown in Table IITTI 

For Si the volume changes AV for an assumed 
first-order /3-tin— > Imma transition from LDA and from 
GGA and are very small compared to the changes for 
Jmmo->sh. Thus we conclude that the phase transition 
/3-tin— > Imma in Si most probably is of second order. 

For Ge the volume changes for a first-order j3- 
tin—tlmma and Imma— >sh transition are within the same 
range. The value of pt = 792 kbar is higher than the ones 
from the estimated range (750-763 kbar) of the previous 
section, and overestimating p t results in overestimating 
AV. A reason for the large value of pt can be found in the 



TABLE III: Transition pressure (pt), order (ord), and volume 
change (AV) for Ge from Fig. E] and similarly for Si. 
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fact, that the last reliable data points are too far away 
from the region where the phase transition is expected to 
occur. Thus, it is not possible to decide whether the tran- 
sition /3-tin— > Imma is of first or of second order. Because 
of the similarities of Si and Ge a second-order transition 
seems more probable. If this phase transition is a discon- 
tinuous one, it is at least weakly of first order. 

Using the enthalpy method the transition pressure for 
cd^/3-tin in Si is calculated to be larger than the one for 
the /3-tin— > Imma transition. These two phase transitions 
take place in a small pressure range and therefore numer- 
ical errors in the determination of the transition pressure 
may lead to a wrong sequence. See the discussion at the 
end of Section IIIII 



VI. ENERGY SURFACES AND ENTHALPY 
BARRIERS 

In this section, we present the results for the energy 
surfaces for Ge (similar calculations have been done for 
Si) as well as for the enthalpy barriers for the cd— > /3-tin 
and Imma— >sh transitions. 

We have obtained the energy surface by calculating the 
energy corresponding to given values of the parameters 
V, b/a, and c/a and a relaxed internal parameter A. On 
the energy surface the pressure is non-hydrostatic except 
one or two lines or a few points. Next, we consider the 
pressure defined in Eq. @ in the non-hydrostatic case 

Po = ~ (Px +Py +Pz) (5) 

where — p x , —p y , and —p z are the components of the 
diagonal stress tensor, Eq.Q. 

For every first-order phase transition, there exists an 
energy barrier which has to be overcome on the path 
from the one phase to the other. For most reactions, 
such as, e.g., adsorption processes, it is a barrier of the 
total energy E which lies between two local minima in 
the corresponding energy space. For this kind of transi- 
tion, the reaction path and the energy barrier between 
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the two phases can be found by using, e.g., the nudged- 
elastic-band method^iH by which the lowest energy 
path between two local minima is detected. However, 
this method requires the existence of well-defined local 
energy minima corresponding to the two phases. 

If the transition, for instance the one from the cd to 
the /3-tin phase in Ge, would occur without any influence 
of pressure, the corresponding energy space would be the 
energy surface E(V, c/a) drawn in Fig. [7| The reaction 
coordinate here is taken as c/a. As expected, there are 
two local minima and a saddle point between them. By 
symmetry, the condition p x = p y is valid all over the en- 
ergy surface. The energy minima and the saddle point 
lie on the contour line p x — p z , which indicates the hy- 
drostatic condition. The energy barrier is then defined 
as the energy difference between the saddle point and the 
starting minimum. 
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FIG. 7: Contour plot of the total energy E(V, c/a) and of the 
average pressure po(V, c/a) (see Eq. for Ge. The interval 
of the contour lines is 50 meV for the energy and 20 kbar 
for the pressure surfaces. The conditions p x = p y , b/a — 1, 
and A = 0.25 are fulfilled within the whole area. The black 
dots mark the equilibrium positions of the cd (c/a = v2) and 
the /3-tin phase (c/a = 0.55). The dotted line markes the 
parameters under hydrostatic condition. 

However, in our case the phase transition is driven by 
the pressure, and therefore the associated thermodynam- 
ical potential is the Gibbs free energy or the enthalpy, if 
temperature effects are neglected. The enthalpy surface 
corresponding to the cd and /3-tin phases for Ge is shown 
in Fig.|Hl The enthalpy for each set of parameters (V, c/a) 
is calculated as H(V,c/a) — E(V,c/a) + V po(V, c/a). In 
the range displayed in Fig. [S] the enthalpy surface has no 
minimum. Therefore, the usual path methods are not 
applicable here. However, the enthalpy barrier can be 
calculated in a different way as follows: The phase tran- 
sition between the A and B phases occurs at the pressure 
p t at which H A (p t ) = H B (p t ). Hence, we have to look 
for the points along the isobars in Fig. |H1 which have the 
same value of the enthalpy. 

In order to visualise the enthalpy barrier, we draw in 
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FIG. 8: Same as Fig.[7|but for the enthalpy H(V,c/a) of Ge. 
The interval between the contour lines for H is 500 meV. 



Fig. [5] the enthalpy at constant average pressure (with 
Px = Py) as a function of the reaction coordinate c/a. 
For convenience the enthalpy is reduced by the starting 
point with c/a = \/2, 



AH(c/a) = H(c/a) - H(V2) 



(6) 



which determines here the cd structure. As can be seen 
in Fig. [5] there exists just one enthalpy line at constant 
pressure, which touches the zero line twice according to 
H A (p) — H B (p) . The corresponding pressure is the tran- 
sition pressure pt and the height of the maximum between 
these minima is the height of the enthalpy barrier for 
this phase transition. Note that on the isobars at most 
three points refer to hydrostatic pressure, while the other 
points corresponds just to the condition p x = p y . 
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FIG. 9: Enthalpy difference at constant average pressure vs. 
the reaction coordinate c/a, see Eq. (|HJ for Ge. The solid line 
is for the cd— >/3-tin transition pressure. The pressure interval 
between consecutive lines is 10 kbar. 

Similar pictures can be drawn for the cd— > /3-tin tran- 
sition of silicon. The height of enthalpy barrier is found 
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to be 423 meV for Ge within GGA and 515 and 508 meV 
for Si within GGA and LDA, respectively. In agreement 
with previous results?^ the barrier for Ge is smaller than 
for Si. 

For the cd— >/?-tin phase transition the determination 
of the enthalpy barrier is easy because only the volume 
and the c/a ratio can vary, and the latter can be used 
as the reaction coordinate. For the /3-tin— >Imma— >sh 
transition all the parameters c/a, b/a, b/c, and A are 
changing (but only two of the ratios are independent). 
We have calculated the total energy for each set of pa- 
rameter (V,c/a,b/a), the internal parameter A has been 
relaxed for each set. Selected energy contour plots at 
constant volume are presented in Fig. ^] In these plots, 
the places where two of the diagonal components of the 
stress tensor are equal are explicitly shown by the dot- 
ted, dashed, and dash-dotted lines. Obviously, hydro- 
static condition is indicated by the points where all of 
these lines cross. Different hydrostatic pressures at one 
volume (like at V — 24 A 3 ) can lead to a multivalued 
pressure- volume relation. 
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FIG. 10: Contour plots of the total energy E(b/a, c/a) for 
selected volumes noted in the insets for Ge. The interval of 
the contour lines is 1 meV. The black dots mark the equilib- 
rium positions. The hydrostatic condition is produced by the 
crossing of all the lines p x — p y , p y = p z , and p z — p x . 

In order to calculate the enthalpy barrier in analogy 
to the cd— >/3-tin transition, we present a contour plot 
of the enthalpy as a function of the volume V and of 
the b/a ratio in Fig. ^] Actually, the enthalpy depends 
also on the other variables c/a and A. The contour plot 
of Fig. ^] has been obtained using the following proce- 
dure: For each value of the set (V,b/a), we have chosen 
c/a so as to fulfill the condition p z = p x (see Fig. 1100 . 
which is suggested because at the end we will look for 
the points which fulfill the hydrostatic condition. The 
value of A has been taken so as to minimise the en- 
ergy for a given choice of the other parameters. The 
two-dimensional presentation of the results as in Fig. 1111 
corresponds to the choice of b/a as the reaction coordi- 



nate. Completely equivalent results can be obtained by 
choosing c/a or A instead. The hydrostatic condition in 
Fig. ^] is then realised by p y = p z . 
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FIG. 11: Contour plot of the enthalpy H(V,b/a) and of the 
pressure po(V,b/a) for Ge. The interval of the contour lines 
is 300 meV for the enthalpy and 20 kbar for the pressure 
increasing from right to left. The condition p x — p z is ful- 
filled within the whole area. The dash-dotted line reproduce 
the hydrostatic condition. The inset shows the schematical 
behavior between 23.5 and 24 A 3 , see text. 



The comparison with our previous results shows that 
most of the data shown in Fig. [5] lie on the hydrostatic 
line. Numerical instabilities lead to deviations (generally 
of less than 1 kbar) between the relaxed points and the 
line p y — p z . Nevertheless, the relaxed points show the 
same behaviour (to a lesser extent) as the hydrostatic 
line, which is shown schematically in the inset of Fig. 1111 

The schematical behaviour of the hydrostatic curve be- 
tween approximately 23.5 and 24 A 3 as shown in the inset 
of Fig.^Jcan be used to discuss the hysteresis effect qual- 
itatively. If one follows the hydrostatic line p y — p z with 
increasing pressure p (decreasing volume) starting from 
A, one arrives at the point C where the direction of the 
line is changing. From C to D, the pressure decreases 
with a simultanous decrease of the volume which is an 
instable situation. Thus, at C the phase transition must 
occur through a jump from C to E, which is at the same 
pressure. From E increasing pressure leads to F. Follow- 
ing the path in the opposite direction starting from F, 
i.e., by decreasing the pressure (volume increase), the di- 
rection of the curve changes at D, which results in a jump 
to B. As a consequence of that it is possible to have a 
higher transition pressure with increasing than with de- 
creasing pressure. Therefore, an alternative estimate for 
the limit of the hysteresis effect can be given. For Si 
within GGA (LDA), we obtaine a transition pressure of 
198 (133) kbar with increasing and of 169 (109) kbar 
with decreasing pressure. For Ge within GGA, the cor- 
responding values are 1010 and 965 kbar. 

A behaviour similar to the one depicted in the inset of 
Fig. ^2 can be found also for the cd— >/?-tin phase tran- 
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sition (see Fig. [HJ. There, the low-pressure edge corre- 
sponding to the point C of the inset of Fig. is at a 
negative pressure, which would lead for /3-tin^cd to a 
very low or even negative transition pressure, the latter 
case being connected to the irreversibility of the phase 
transition. For the cd— >/?-tin transition the volume and 
the enthalpy as a function of hydrostatic pressure is pre- 
sented in Fig. The values have been extracted from 
Fig. [S]along the hydrostatic line p x = p z (= p y ). The ideal 
cd structure (c/a = y/2) has been reached within an er- 
ror of 1% and therefore cd is noted as CD in analogy 
to the difference between sh and SH (see Section . 
In order to discriminate the enthalpy curves of Fig. IT51 
we have subtracted a linear background to arrive at the 
curves AH. The solid and the dashed line mark the CD 
and the BCT structure, respectively. The dotted line in- 
dicate the values along the part of the hydrostatic line, 
which connects the CD and the BCT structure. The des- 
ignated points correspond to the inset of Fig. 1111 
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FIG. 12: Volume, enthalpy and reduced enthalpy vs. hydro- 
static pressure for Ge. The points noted in the figures corre- 
spond to the points of the inset of Fig. 1111 T marks the point 
of the phase transition between the cd and the /3-tin phase, 
see text. 

The condition of local stability is fulfilled along all 
curves shown in Fig. 1121 Thus, the volume vs. pressure 
curves are strictly monotonically decreasing and the en- 
thalpy vs. pressure curves are convex. The condition of 
global stability requires the enthalpy to be minimal, and 
so the curves from A to T and from T to F mark the 
globally stable regime. The point of the phase transi- 
tion in analogy to the Maxwell construction is noted as 



T. The pressure range of the coexistence regime is de- 
termined by the overlap of the globally stable and the 
instable regions, which is here from B (or D) to E (or C). 

A similar picture could be drawn for the Jmma-»sh 
transition, but with less resolution. A reliable identifica- 
tion of the structures within the H(p) curve analoguos 
to Fig. 1121 requires a finer grid of data points for the in- 
terpolation and/or a higher convergence. 

The enthalpy barrier for the cd— >/3-tin phase transition 
can be determined with a procedure like illustrated in 
Fig. In Figs. EI] and [21 we show the enthalpy at 
given pressures calculated within GGA as a function of 
the reaction coordinate b/a, for Si and Ge, respectively. 



x 

< 




p =50 kbar 

— p = 1 94 kbar 

- - 50<p <250 kbar 



0.90 0.92 



0.94 0.96 
b/a 



0.9E 



1 .00 



FIG. 13: Enthalpy difference at constant average pressure 
vs. reaction coordinate b/a, see Eq. J7J for Si. The pressure 
interval between the dashed lines is 20 kbar. The solid line 
marks the transition pressure found here. 
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FIG. 14: Same as Fig. ED but for Ge. The pressure interval 
between the dashed lines is 40 kbar. The solid curve marks the 
line corresponding to the transition pressure for the /mma-tsh 
transition obtained in Section Ivl 

The enthalpy in Figs.ll3landll4lis reduced by the value 
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at b/a — 1 (corresponding to the /3-tin phase), 

AH(b/a) = H{b/a) - H(l) . (7) 

The dotted line in Fig. ll3l corrcsponds to a pressure where 
the enthalpy has a local minimum for the /3-tin phase. At 
a pressure of about 130 kbar we find a continuous transi- 
tion to the Imma phase. Above this pressure, a minimum 
appears at b/a < 1. By further increasing po a second 
minimum apprears at b/a ~ 0.92, corresponding to the sh 
phase. At the pressure of 194 kbar the two minima have 
the same enthalpy, therefore a transition can occur. The 
enthalpy barrier between the two phases is surprisingly 
small, namely 1.2 meV. At higher pressure the sh phase 
is stable. For the LDA calculation a transition pressure 
of 133 kbar and an enthalpy barrier of 1.3 meV are ob- 
tained. If one considers the enthalpy curves of Fig. H~3l 
at the transition pressures determined in the Sections IIVI 
and one can give a further estimate of the enthalpy 
barriers, and, the values obtained following this proce- 
dure do not differ appreciably from the ones mentioned 
above. 

As can be seen in Fig. for Ge, a minimum of the 
enthalpy can be identified corresponding to a pure Imma 
phase, slightly below 820 kbar. However, within our res- 
olution (about 0.1 meV) the pressure curves have only 
one minimum. Therefore, an enthalpy barrier between 
the Imma and the sh phase can not be identified. In- 
deed, this transition is of first order. Hence, the existing 
barrier must be smaller than 0.1 meV. 

All the enthalpy barriers which we are able to obtain 
for the lmma->sh transition are much smaller than the 
thermal energy at room temperature (25 meV). Thus, 
experiments performed at room temperature should not 
be able to distinguish between a first- or a second-order 
phase transition. 

VII. DISCUSSION OF THE RESULTS FOR THE 
PHASE TRANSITIONS 

The transition pressures for the cd^/3-tin^/mma— >sh 
phase transitions have been calculated with different 
methods in the previous sections. In Sections IIVI and IVT1 
an estimate of the hysteresis range was given by consider- 
ing the relaxed values (R) and the hydrostatic line (HL) 
respectively. In Section [V] the transition pressure was 
determined via the enthalpy vs. pressure curves (EC). In 
Table IIVI an overview of all our results for the transition 
pressure is given in comparison with the experimental 
results. 

For Si the value of the pressure for the cd— > /3-tin tran- 
sition calculated with GGA is within the range of the 
experimental values whereas the LDA calculation under- 
estimates it. We could not reproduce the experimentally 
observed hysteresis^ for the /3-tin— ► Imma transition in Si. 
For a first-order phase transitions, such as Imma— >sh, a 
hysteresis could occur due to the behaviour of the hy- 
drostatic line as in Fig. 1111 Since our results support 



TABLE IV: Summary of the estimated transition pressures 
(pt) between the cd, /3-tin, Imma, and sh phase for increas- 
ing (|) and decreasing (1) pressure, determined by the en- 
thalpy curves (EC), the relaxation (R), and the hydrostatic 
line (HL). 
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a second-order phase transition for /3-tin— ► Imma rather 
than a first-order one, no hysteresis should be found in 
this case. The transition pressures mentioned in Ref. 
are obtained either from the first appearance of a new 
phase or from the complete transition of the sample from 
one phase to the other. Thus we estimate the experi- 
mental range in which the /3-tin— > Imma and Imma->sh 
transitions occur in analogy to the procedure explained 
in Section IIVI Unfortunately, the calculated transition 
pressures for /3-tin— ► Imma lie below the ones for the first 
cd^/3-tin phase transition and underestimate the exper- 
imental ones. 

The next phase transition to be considered is 
Imma— >sh. Here a hysteresis was found within our calcu- 
lations as well as in the experiment. The results for the 
transition pressures within LDA underestimate again the 
experimental ones, whereas the ones within GGA over- 
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estimate them. Not included in Table IIVI are the exper- 
imental pressures for mixed (50:50) Imma-sh samples, 
which are 154 kbar for pressure increase and 159 kbar 
for pressure decrease^. This is the only case in which the 
transition to the (mixed) Imma phase occur at a lower 
pressure for an increasing than for a decreasing external 
pressure. All other experimental points show the same 
kind of hysteresis which is reproduced in our results. The 
difference of the pressures, where the Imma phase was 
measured with pressure increase and decrease is 5 kbar 
for the mixed samples and 11-12 kbar for the pure ones. 
From the experimental transition pressures the hysteresis 
is found to be between 3 and 14 kbar. Within our results 
the largest possible hysteresis is found to be around 29 
(24) kbar and the smallest around 1-12 (2-7) kbar for 
GGA (LDA), which is in agreement with the experimen- 
tal values. The transition pressures determined by the 
crossing of the enthalpy curves are found to be within 
the hysteresis range calculated here. 

For Ge the transition pressure for the cd^/3-tin tran- 
sition is slightly smaller than the experimental one but 
the results for the /3-tin— >-/mma transition match much 
better. In the experiment^* just the /3-tin— >Imma tran- 
sition was examined, because the transition pressure for 
the Imma— >sh transition was higher than the maximum 
accessible pressure of 810 kbar. In comparison to the ex- 
perimental results of Si it is assumed that the Imma phase 
for Ge is stable from 750 to 850 kbar, which is compatible 
with the observation of a sh phase at 900 kbari^i From 
this point of view, the transition pressure calculated here 
is slightly too high. 

There are many reasons for the discrepancy between 
our results and the experimental ones. On one hand, the 
experimental conditions like pressure profiles or measure- 
ment times are not perfectly controlled, whereas nonhy- 
drostaic conditions affect the transition pressures; on the 
other hand, temperatue effects, which are neglected in 
our calculation, could be of some relevance. In fact, for 
Si the theoretical transition pressure decreases with in- 
creasing temperature^ 



VIII. CONCLUSIONS 

We have presented a first-principles investigation of the 
pressure-induced phase transitions cd — > /3-tin — > Imma 
— > sh in Si and Ge. Numerical accuracy has required 
these calculations to be performed within the same unit 
cell and the pressure to be calculated from the stress 
tensor instead of using any equation of state. We have 
determined the phase transitions by investigating the re- 
laxed parameters, the enthalpy vs. pressure curves, and 
the hydrostatic condition. From the pressure dependence 
of the relaxed structure parameters the hysteresis effect 
does not seem to have been verified in a theoretical inves- 
tigation before. An examination of the hydrostatic con- 
dition provides a thermodynamical explanation for the 
existence of this hysteresis. Further on, the orders of the 



phase transitions have been investigated with the helps of 
the curves of the enthalpy vs. pressure and of the deriva- 
tive of the enthalpy vs. pressure. For Si the transitions 
Imma— >sh and /3-tin— > Imma are found to be of first and 
of second order, respectively. For Ge the results support 
a second-order /3-tin— > Imma transition more than a first- 
order one, but in this case the order cannot be definitely 
determined. A calculation of the enthalpy barriers be- 
tween the phases show a vanishing barrier between the 
/3-tin and the Imma phase for Si, which also supports the 
second-order /3-tin— > Imma transition, whereas no barrier 
could be found for all /3-tin— > Imma— >sh transitions. Fi- 
nally, our alternative results are fully consistent with each 
other and show an acceptable agreement with the avail- 
able experimental data. 
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APPENDIX A: FERMI ENERGY AND BAND 
STRUCTURE 

We examined the volume dependence of the Fermi en- 
ergy of BCT, BCO, and SH. In Fig^we show the Fermi 
energy reduced by the BCO values for Ge, 
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FIG. 15: Fermi energies of the involved structures for Ge 
reduced by the Fermi energy of BCO, taken as zero. 

In the range of stability of the /3-tin phase the Fermi 
energy of BCO is identical with the one of BCT. Simi- 
larly, the Fermi energy of BCO and SH are identical in 
the range of stability of the sh phase. Therefore, the 
regime where the Fermi energy of BCO is neither equal 
to the one of BCT nor to the one of SH is assumed here 
to be the range of stability of the pure Imma phase. 
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We analysed the volume dependence of the band struc- 
ture for the BCT, BCO, and SH structure, too. The re- 
sults for the lowest four bands at the T-point for Ge are 
displayed in Fig. 

The energy values Ei are reduced in such a way that 
the Fermi energy at any volume is set equal to zero, 

AE(V) = Ei(V) — E F (V) . (A2) 

In most of the volume range presented in Fig. 1161 the 
band structure presents a linear dependence on volume. 
However, the energy gradients differ for each structure. 
Within the transition region between the Imma and the 
sh phase one observes a small discontinuity (here appar- 
ent at V = 23.64A 3 ) which is attributed to the first- 
order phase transition. When starting from the /3-tin 
phase, a deviation of the BCO band energies from the 
BCT energies can be observed, which increases contin- 
uously. Thereby, the degeneracy of the BCT bands is 
lifted, which indicates the appearance of the less symmet- 
rical Imma phase. The behavior of the band structure for 
Si is very similar to the one of Ge. 
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FIG. 16: Eigenenergies at the V point, reduced by the Fermi 
energy for the BCO, BCT and SH structure for Ge. 
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